/****************************************************************
 * Fulguro
 * Copyright (C) 2004 Christophe Clienti
 *
 * This library is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Lesser General Public
 * License as published by the Free Software Foundation; either
 * version 2 of the License, or (at your option) any later version.
 *
 * This library is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
 * Lesser General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public
 * License along with this library; if not, write to the
 * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
 * Boston, MA 02111-1307, USA.
 ***************************************************************/


#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <flgrCoreVector.h>
#include <flgrCoreArithVector.h>
#include <flgrCoreDispatch.h>
#include <flgrCoreDataIO.h>
#include "flgrMatrixTransform.h"
#include "flgrMatrixTransformSSE.h"
#include "flgrMatrixTransformSSE2.h"

////////////////////////////////////////////////////////////////////////////
// macro for different type of multdata_dtype
////////////////////////////////////////////////////////////////////////////
#define flgr2d_matrix_matrixprod(dtype)					\
  int i,j,k;								\
  int n,n2,m,m2;							\
  FLGR_Vector *vdat_a = flgr_vector_create(dat_a->spp,dat_a->type);	\
  FLGR_Vector *vdat_b = flgr_vector_create(dat_a->spp,dat_a->type);	\
  FLGR_Vector *vsum = flgr_vector_create(dat_a->spp,dat_a->type);	\
									\
  FLGR_DECLARE_FUNCTION;						\
									\
  m=dat_a->size_y;							\
  n=dat_b->size_x;							\
  m2=dat_b->size_y;							\
  n2=dat_a->size_x;							\
									\
  for(i=0 ; i<m ; i++) {						\
    for(j=0 ; j<n ; j++) {						\
      flgr_vector_clear(vsum);						\
      for(k=0 ; k<m2; k++) {						\
	flgr2d_get_data_vector_##dtype(dat_a,i,k,vdat_a);		\
	flgr2d_get_data_vector_##dtype(dat_b,k,j,vdat_b);		\
	flgr_vector_mult(vdat_a, vdat_a, vdat_b);			\
	flgr_vector_add(vsum,vsum,vdat_a);				\
      }									\
      flgr2d_set_data_vector_##dtype(dat_res,i,j,vsum);			\
    }									\
  }									\
									\
  flgr_vector_destroy(vdat_a);						\
  flgr_vector_destroy(vdat_b);						\
  flgr_vector_destroy(vsum);						\
									\
  return

void flgr2d_matrix_mult_fgBIT(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr_no_define_type_function(FLGR_BIT);
}

void flgr2d_matrix_mult_fgUINT8(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgUINT8);
}

void flgr2d_matrix_mult_fgUINT16(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgUINT16);
}

void flgr2d_matrix_mult_fgUINT32(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgUINT32);
}

void flgr2d_matrix_mult_fgINT8(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgINT8);
}

void flgr2d_matrix_mult_fgINT16(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgINT16);
}

void flgr2d_matrix_mult_fgINT32(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgINT32);
}

void flgr2d_matrix_mult_fgFLOAT32(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgFLOAT32);
}

void flgr2d_matrix_mult_fgFLOAT64(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  flgr2d_matrix_matrixprod(fgFLOAT64);
}

FLGR_Ret flgr2d_matrix_mult(FLGR_Data2D *dat_res, FLGR_Data2D *dat_a, FLGR_Data2D *dat_b) {
  int n,m;
  FLGR_DECLARE_FUNCTION;
  
  if((dat_res==NULL) || (dat_a==NULL) || (dat_b==NULL)){
    POST_ERROR("Null objects!\n");
    return FLGR_RET_NULL_OBJECT;
  }

  m=dat_a->size_y;
  n=dat_b->size_x;

 
  if((m!=n)) return -1;
  if((dat_res->size_y!=m) || (dat_res->size_x!=n)) return FLGR_RET_SIZE_ERROR;
  if(dat_a->type != dat_b->type) return FLGR_RET_TYPE_DIFFERENT;
  if(dat_res->type != dat_b->type) return FLGR_RET_TYPE_DIFFERENT;

  FLGR_DISPATCH_PROCEDURE(dat_a->type,flgr2d_matrix_mult,dat_res, dat_a, dat_b);

}
